rm(list=ls())
library(foreign)
library(locfit)

lgraph = function(X,Y,CONDITION, main="", xlab="X", ylab="Y", ylim=NULL, deg=2, ...){
  fit <- locfit(Y[CONDITION]~lp(X[CONDITION], deg=deg, ...)); crit(fit) <- crit(fit,cov=0.95); 
  plot(fit,band="local", main=main, xlab=xlab, ylab=ylab, ylim=ylim)
}
LGBT=read.dta("~/Dropbox/LGBT/data/output/LGBT_long.dta", convert.factors=FALSE) 
LGBT=as.data.frame(LGBT)

LGBTM=aggregate(cbind(lgbt_articles, mlgbt_articles) ~ year, LGBT, mean) 
# keep obs over 2002 given the spareseness of data in 2000-2
LGBTM = LGBTM[LGBTM$year>2002,]

#######################################################################################################
# Figure 1: LGBT Saliency 
#######################################################################################################
pdf(file="~/Dropbox/LGBT/draft/figures/LGBTSaliencyOT.pdf", width=12, height=5) 
deg=3
par(mfrow=c(1,2))
lgraph(LGBTM$year, LGBTM$lgbt_articles, TRUE, main="Highest Circulation Newspaper", xlab="", ylab="Mean Number of LGBT Articles", deg=deg, ylim=c(5, 35))
lgraph(LGBTM$year, LGBTM$mlgbt_articles, TRUE, main="Mean Available Newspapers", xlab="", ylab="Mean Number of LGBT Articles", deg=deg, ylim=c(5, 35))
dev.off()

#######################################################################################################
# Figure 2: LGBT issue saliency by country 
#######################################################################################################

require(ggplot2)
LGBT=LGBT[LGBT$year>2002,]
LGBT=LGBT[LGBT$countryname!="Malawi",]
pdf(file="~/Dropbox/LGBT/draft/figures/LGBTSaliencybyCountry.pdf", width=12, height=12) 
f1 = ggplot(LGBT, aes(x=year, y=mlgbt_articles))+geom_point(aes(color=countryname))
f1 = f1 + facet_wrap(~countryname)+ labs(title="Country level trajectries (2003-2013)", y="Annual LGBT articles (mean available newspaprs)", x="Year") +scale_x_continuous(breaks=seq(from=2003, to=2013, by=2)) 
f1 = f1 + theme(axis.text.x=element_text(angle=90, vjust=1, hjust=0)) 
print(f1)
dev.off()


#######################################################################################################
# Figure 4: LGBT saliency and key explanatory variables
#######################################################################################################

# Renewalist Pop Share 
medrs=median(LGBT$renwshare_pop, na.rm=T)
LGBT$renwshareH=as.numeric(with(LGBT, LGBT$renwshare_pop>medrs))

RENWY=aggregate(cbind(lgbt_articles, mlgbt_articles) ~ year+ renwshareH, LGBT, mean, na.rm=T) 
RENWY$renwshareHN=RENWY$renwshareH
RENWY$renwshareHN=ifelse(RENWY$renwshareHN==0, 1,16)
RENWY=RENWY[RENWY$year>2002,]
print(RENWY)

# Polity IV  
medpol=median(LGBT$polity2, na.rm=T)
LGBT$democH=as.numeric(with(LGBT, LGBT$polity2>medpol))

DEMY=aggregate(cbind(lgbt_articles, mlgbt_articles) ~ year+ democH, LGBT, mean, na.rm=T) 
DEMY$democHN=DEMY$democH
DEMY$democHN=ifelse(DEMY$democHN==0, 1,16)
DEMY=DEMY[DEMY$year>2002,]

DEMY=DEMY[DEMY$year>2002 & DEMY$year<2012,]
print(DEMY)

# Opposition Seat Share (Margin of Minority)
mePolitcomp=median(LGBT$Politcomp, na.rm=T)
LGBT$PolitcompH=as.numeric(with(LGBT, LGBT$Politcomp>=mePolitcomp))
table(LGBT$PolitcompH)

PolitcompH=aggregate(cbind(lgbt_articles, mlgbt_articles) ~ year+ PolitcompH, LGBT, mean, na.rm=T) 
PolitcompH$PolitcompHN=PolitcompH$PolitcompH
PolitcompH$PolitcompHN=ifelse(PolitcompH$PolitcompHN==0, 1,16)
PolitcompH=PolitcompH[PolitcompH$year>2002,]
print(PolitcompH)

# State Legitimacy
medlegit=median(LGBT$legit, na.rm=T)
LGBT$LegitH=as.numeric(with(LGBT, LGBT$legit>=medlegit))
LEGITY=aggregate(cbind(lgbt_articles, mlgbt_articles) ~ year+ LegitH, LGBT, mean, na.rm=T) 
LEGITY$LegitHN=LEGITY$LegitH
LEGITY$LegitHN=ifelse(LEGITY$LegitHN==0, 1,16)
LEGITY=LEGITY[LEGITY$year>2002,]
print(LEGITY)

pdf(file="~/Dropbox/LGBT/draft/figures/LGBTSaliencyOverTimePanelB.pdf", width=10, height=9) 
deg=3
par(mfrow=c(2,2), mar = c(4, 3, 2, 1))
plot(RENWY$lgbt_articles~RENWY$year, xlab="", ylab="LGBT articles", 
     pch=RENWY$renwshareHN, ylim = c(0, 50), axes=FALSE, cex=2)
axis(2, at = c(0, 25, 50), cex.axis=2)
axis(1, at = c(2003, 2007, 2011), cex.axis=2)
title(main="Renewlaist pop share", cex.main=2)
H.Loess <- loess(lgbt_articles[renwshareH==1]~year[renwshareH==1], data = RENWY)
FitH <- fitted(H.Loess)
lines(RENWY$year[RENWY$renwshareH==1], FitH, col="dark blue", lt=2)
L.Loess <- loess(lgbt_articles[renwshareH==0]~year[renwshareH==0], data = RENWY)
FitL <- fitted(L.Loess)
lines(RENWY$year[RENWY$renwshareH==0], FitL, col="dark red")

plot(PolitcompH$lgbt_articles~PolitcompH$year, xlab="", ylab="LGBT articles", 
     pch=PolitcompH$PolitcompHN, ylim = c(0, 50), axes=FALSE, cex=2)
axis(2, at = c(0, 25, 50), cex.axis=2)
axis(1, at = c(2003, 2007, 2011), cex.axis=2)
title(main="Opposition Seat Share", cex.main=2)
H.Loess <- loess(lgbt_articles[PolitcompH==1]~year[PolitcompH==1], data = PolitcompH)
FitH <- fitted(H.Loess)
lines(PolitcompH$year[PolitcompH$PolitcompH==1], FitH, col="dark blue", lt=2)
L.Loess <- loess(lgbt_articles[PolitcompH==0]~year[PolitcompH==0], data = PolitcompH)
FitL <- fitted(L.Loess)
lines(PolitcompH$year[PolitcompH$PolitcompH==0], FitL, col="dark red")

plot(DEMY$lgbt_articles~DEMY$year, xlab="", ylab="LGBT articles", 
     pch=DEMY$democHN,  ylim = c(0, 50), axes=FALSE, cex=2)
axis(2, at = c(0, 25, 50),  cex.axis=2)
axis(1, at = c(2003, 2007, 2011),  cex.axis=2)
title(main="Polity IV", cex.main=2)
H.Loess <- loess(lgbt_articles[democH==1]~year[democH==1], data = DEMY)
FitH <- fitted(H.Loess)
lines(DEMY$year[DEMY$democHN==1], FitH, col="dark blue", lt=2)
L.Loess <- loess(lgbt_articles[democH==0]~year[democH==0], data = DEMY)
FitL <- fitted(L.Loess)
lines(DEMY$year[DEMY$democH==0], FitL, col="dark red")

plot(LEGITY$lgbt_articles~LEGITY$year, xlab="", ylab="LGBT articles", 
     pch=LEGITY$LegitHN,  ylim = c(0, 50), axes=FALSE, cex=2)
axis(2, at = c(0, 25, 50),  cex.axis=2)
axis(1, at = c(2003, 2007, 2011),  cex.axis=2)
title(main="State Legitimacy", cex.main=2)
H.Loess <- loess(lgbt_articles[LegitH==1]~year[LegitH==1], data = LEGITY)
FitH <- fitted(H.Loess)
lines(LEGITY$year[LEGITY$LegitH==1], FitH, col="dark blue", lt=2)
L.Loess <- loess(lgbt_articles[LegitH==0]~year[LegitH==0], data = LEGITY)
FitL <- fitted(L.Loess)
lines(LEGITY$year[LEGITY$LegitH==0], FitL, col="dark red", lt=1)
dev.off()


#######################################################################################################
# Appendix Figure 2: Legislative Events 
#######################################################################################################

library(plyr)
myvars <- c("legalevent", "mlgbt_articles")
legal <- LGBT[myvars]

legal=legal[!is.na(legal$legalevent), ]
legal=legal[!is.na(legal$mlgbt_articles), ]
mm <- ddply(legal, "legalevent", summarise, mlgbt = mean(mlgbt_articles))
mm$legalevent <- factor(mm$legalevent,
                        levels = c(0,1,2),
                        labels = c("No event", "Law debated", "Law passed"))

pdf(file="~/Dropbox/LGBT/draft/figures/LegalEvents.pdf", width=12, height=12) 
f2= ggplot(mm, aes(x = factor(legalevent), y = mlgbt)) + geom_bar(stat = "identity")
f2 = f2 + labs(title="", y="Mean Number of LGBT Articles", x="") 
f2 = f2 + theme(axis.text.x = element_text(colour="grey20",size=25,hjust=.5,vjust=.5,face="plain"), 
                axis.text.y = element_text(colour="grey20",size=18,hjust=1,vjust=0,face="plain"),
                axis.title.y = element_text(colour="grey20",size=22,hjust=.5,vjust=.5,face="plain"))
print(f2)
dev.off()